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In this study we model early times dynamics of relativistic heavy ion collisions by an initial color 
electric field which then decays to a plasma by the Schwinger mechanism, coupling the dynamical 
evolution of the initial color field to the dynamics of the many particles system produced by the 
decay. The latter is described by relativistic kinetic theory in which we fix the ratio rj/s rather 
than insisting on specific microscopic processes, and the backreaction on the color field is taken into 
account by solving self-consistently the kinetic and the field equations. We study isotropization and 
thermalization of the system produced by the field decay for a static box and for a 1 + ID expanding 
geometry. We find that regardless of the viscosity of the produced plasma, the initial color electric 
field decays within 1 fm/c; however in the case 77 /s is large, oscillations of the field are effective along 
all the entire time evolution of the system, which affect the late times evolution of the ratio between 
longitudinal and transverse pressure. In case of small 77 /s ( 77 /s < 0.3) we find nsotropization ~ 0.8 
fm/c and Tthermaiization ~ 1 fm/c in agreement with the common lore of hydrodynamics. Moreover 
we have investigated the effect of turning from the relaxation time approximation to the Chapman- 
Enskog one: we find that this improvement affects mainly the early times evolution of the physical 
quantities, the effect being milder in the late times evolution. 
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I. INTRODUCTION 


The understanding of early times dynamics is one of 
the most interesting and compelling problems of heavy 
ion collisions at ultrarelativistic energy. According to 
the standard picture of such collision processes, imme¬ 
diately after the collision an ensemble of strong longi¬ 
tudinal color electric and color m^netic coherent fields 
known as Glasma is produced [iH^- Inhomogeneities in 
the Glasma fields are amplified by plasma instabilities, 
which eventually lead to loss of coherence of the color 
fields and to a locally isotropized particle plasma [ffil. 
Recently the early isotropization problem has been stud¬ 
ied also by means of AdS/CFT methods (2214^ . 


Besides plasma instabilities, a mechanism responsi¬ 
ble for the initial field decays is the one introduced by 
Schwinger in the context of Quantum Electrodynamics 
[ttI [ 2 ^ , known as the Schwinger effect which consists in 
a vacuum instability towards the creation of particle pairs 
by a strong electric field, and it is related to the existence 
of an imaginary part in the quantum effective action of 
a pure electric field, see for a review. The prob¬ 
lem of pair formation in strong electric fields has been 
considered recently by means of real time lattice simu¬ 
lations (sol - ls^ . Moreover non abelian generalizations of 
the Schwinger production rate has been found for the 
case of both static and time dependent field ■ In 

the context of quark-gluon plasma physics it has been 
considered as a mechanism for the color field decays in 


m. 


The purpose of the present study is to model early 
times dynamics of the system produced in relativistic 


heavy ion collisions by an initial color electric field which 
then decays to a plasma by the Schwinger effect. The 
focus of our study is to compute quantities which serve 
as indicators of thermalization and isotropization of the 
plasma. In particular, we couple the dynamical evolution 
of the color field to the dynamics of the many particles 
system produced by the decay, the latter being described 
by relativistic kinetic theory [6T@- As in the previous 
studies on this subject we assume that the dynamics of 
the color field is abelian, hence it satisfies the classical 
field Maxwell equations, see for a review. Moreover 
in the present approach we neglect for simplicity initial 
longitudinal color magnetic fields which are present in the 
more complete description of the Glasma state. We con¬ 
sider two simple cases: the first one is a system evolving 
in a static box; then we consider a system with a longi¬ 
tudinal expansion which has a greater physical interest 
than the static box because it is closer to the picture of 
the early times dynamics of relativistic heavy ion colli¬ 
sions. 

Besides being the first study in which a Monte Carlo 
method is used to simulate the Schwinger effect in the 
context of early times dynamics of high energy collisions, 
we improve previous studies which mainly rely on Re¬ 
laxation Time Approximation (RTA) [^, or on a 
linearization of the conductive electric current (45l| . by 
avoiding any ansatz both on the electric current and on 
the collision integral in the Boltzmann equation which 
permits to go beyond the RTA. The approach developed 
here, based on a stochastic solution of the relativistic 
Boltzman equation, has in perspective the advantage to 
be easily extended to 3 -I- ID realistic simulations, allow¬ 
ing to study the possible impact of the early dynamics 
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on observables like radial and elliptic flow [^ . 

In both the cases of the static box and of the longitu¬ 
dinal expanding system we have focused on the calcula¬ 
tion of Pl/Pt, where Pl and Pt denote the longitudinal 
and transverse pressures respectively. For the particu¬ 
lar initial field configuration used in this work, namely 
an electric field along the longitudinal direction, one has 
Pl/Pt = — 1 ; however the field decays onto particles 
thanks to the Schwinger effect and particle dynamics re¬ 
moves this initial anisotropy, the efficiency of isotropiza- 
tion being related to the coupling among the particles. 
In our simulations we fix the ratio 77/s where 77 is the 
shear viscosity and s the entropy density, and compute 
locally the corresponding cross section hy means of the 
Chapman-Enskog (CE) approximation . We find 

that increasing 77/s the system is less efficient in removing 
the initial anisotropy in agreement with previous studies, 
although regardless of 77/s the longitudinal pressure be¬ 
comes positive within a small fraction of fm/c, r Ri O.I 
fm/c; moreover for large values of 77/s we find persistent 
oscillations of the color field which affects the evolution 
of physical quantities like PjJPt and the proper energy 
density, in agreement with [SJ, . Even if we go beyond 
the RTA we find that only the early times evolution is af¬ 
fected, leaving the main behaviours very mildly affected. 

The plan of the article is as follows. In Section II 
we review briefly the abelian flux tube model which we 
implement in our simulations. In Section III we review 
relativistic transport theory. In Section IV we present 
our results of the simulations for the case of a static box. 
In Section V we discuss in detail our results for the case 
of a box with a longitudinal expansion. Finally in Section 
VI we draw our conclusions. 


II. ABELIAN FLUX TUBE MODEL 

In this Section we briefly summarize the abelian flux 
tube model (AFTm) [ 39 l - l 5 fil| which we implement in our 
simulations of the initial stage of relativistic heavy ion 
collisions. In the present study we do not insist on im¬ 
plementing the most realistic geometrical condition rel¬ 
evant for heavy ion collisions experiments, in which one 
should take into account several flux tubes in the trans¬ 
verse plane similar to the glasma configuration; rather we 
consider a simpler situation in which there is only one flux 
tube of a given transverse area, and study its dynamical 
evolution by coupling the field equations to relativistic 
kinetic theory for the particle quanta produced by the 
decay of the field itself, leaving to upcoming works the 
study of more realistic initial conditions. 

The main assumptions of the AFTm are: 

• In the initial condition a color electric field is 
present, which is produced by color charges in the 
colliding nuclei; 

• The color electric field decays into particle quanta 
by the Schwinger mechanism; 


• The particle quanta propagate in the medium col¬ 
liding and interacting with the background of the 
color field; 

• Particle creation as well as particle currents affect 
in a self-consistent way the color electric field; 

• The field dynamics is abelian, namely it satisfies 
the abelian Maxwell equations. 

The last assumption is quite strong because initial gauge 
fields are quite large, nevertheless it permits the easiest 
implementation of the coupling of classical field equations 
to relativistic kinetic theory in presence of the Schwinger 
mechanism, hence we will limit ourselves to consider 
abelian dynamics of the classical fields leaving the intro¬ 
duction of non abelian field dynamics to future works. 

We will discuss more about kinetic theory in the next 
Section, therefore in this Section we mainly focus on the 
Schwinger mechanism which is responsible for particle 
production in the early stage of the collision, and on 
the Monte Carlo implementation we use to simulate the 
Schwinger effect. 

In this Section we closely follow at the same time 
adopting a slightly different notation which is more con¬ 
venient for this article. Moreover, since the very begin¬ 
ning we assume the initial color field to be polarized along 
the 3 ’’'^ direction of adjoint color space, meaning that only 
one particular color charge is present in the two colliding 
nuclei [H^. Moreover we assume only one component of 
the electric field is present, namely the one in the longi¬ 
tudinal direction, which we denote by E in the following. 
The latter is a consequence of the fact that in the initial 
condition we assume the field is purely longitudinal, and 
transverse currents are not produced during time evolu¬ 
tion of the system if transverse expansion does not take 
place implying a vanishing transverse field. The number 
of pairs per unit of spacetime and invariant momentum 
space produced by the decay of the electric field by the 
Schwinger effect is (we assume massles quanta through¬ 
out this article) 


dNjc _ dNjc 

dr d^xd^prdpz 


'P-jc{PT)d{pz)PO , 


( 1 ) 


with 




( 2 ) 


the plus (minus) sign corresponding to the creation of a 
boson (fermion-antifermion) pair. In this equation px, 
Pz refer to each of the two particles created by the tun¬ 
neling process; Ejc is the effective force which acts on the 
tunneling pair and it depends on color and flavor; it can 
be written as 


— (sIQjcA'I '^i) ^ (fflQjcA'I CTj) , (3) 

where Uj denotes the string tension depending on the 
kind of flavor considered. Moreover po = \/'Pt Pz cor¬ 
responds to the single particle kinetic energy. 
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The Qjc are color-flavor charges which, in the case of 
quarks, correspond to the eigenvalues of the T3 operator: 

Qji = 2 ’ = “ 2 ’ ^ ^ 

for antiquarks, corresponding to negative values of j, 
the color-flavor charges are just minus the correspond¬ 
ing charges for quarks. Finally for gluons (which in our 
notation correspond to j = 0) the charges are obtained 
by building gluons up as the octet of the 3 0 3 in color 
space; in particular 

Qoi = 1 , Qo2 = 2 J Qo3 ~ ~ 2 ^ 

and Qo4 = —Qoii Q05 = —Q02, Qoe = ~Qo3- The above 
charges are obtained easily [^: for example, the gluon 
we label with subscript 1 can be represented, from the 
point of view of transformations in color space, as a red- 
antigreen state; its T3 charge is then given by difference 
of T3 charges of a red and a green state. We notice that 
we have only six gluons out of eight, corresponding to 
the non-diagonal color generators; the gluon fields corre¬ 
sponding to the diagonal color generators have vanishing 
coupling with the background field, hence they cannot be 
produced by the Schwinger effect. 

As a final remark we would like to observe that, even 
if the model is named abelian, such nomenclature sim¬ 
ply refers to the fact that in the evolution equation for 
the classical field, self-interacion terms coming from non 
vanishing structure constants of the color gauge group 
are neglected [s^. However interactions among the clas¬ 
sical field and gluons are still present in this calcula¬ 
tions, thanks both to the Schwinger effect which produces 
charged gluons, and to conduction currents which affect 
the evolution of the field, see the next Section for more 
details. 


III. RELATIVISTIC TRANSPORT THEORY 
COUPLED TO MAXWELL EQUATIONS 

Our calculation scheme is based on the Relativistic 
Transport Boltzmann equation which, in the presence of 
a gauge field F^’^, can be written as follows: 

+ gQjcF^^''pudl) fjc{x,p) = -h Cjclf] (6) 

where fjc{x,p) is the distribution function for flavour j 
and color c, is the electromagnetic tensor. On the 
right hand side we have the source term dN/dT which 
describes the creation of quarks, antiquarks and gluons 
due to the decay of the color electric field and C[f] which 
represents the collision integral. Considering only 2 —>■ 
2 body elastic scatterings, the collision integral can be 
written as: 

C\f] = f r r r. 

7 2F2(27r)3 2Fr(27r)3 2F'(27r)3 

x\M\‘^S'^{pi+P2-pi'-P2>) , ( 7 ) 


where we omit flavour and color indices for simplicity, M 
is the transition matrix for the elastic process linked to 
the differential cross section through I Afp = 16ns'^da/dt, 
being s the Mandelstam variable. 

In our simulations we solve numerically Eq. m using 
the test particles method and the collision integral is com¬ 
puted using Monte Carlo methods based on the stochas¬ 
tic interpretation of transition amplitude [dOl - l^ . [^ . . 

The evolution of the electric field is given by the 
Maxwell equations: 


dE _ dE _ 


( 8 ) 


where p corresponds to the color charge density, 


P = 


9'^ Qjc j d^Pfjcip) , 

j,c 


( 9 ) 


where j, c denote flavor and color respectively; the sum in 
the above equation runs over quarks, antiquarks and glu¬ 
ons. On the other hand j corresponds to the color electric 
current which is given by the sum of two contributions: in 
fact the Schwinger effect can be described as a dielectric 
breakdown in which dipoles are produced by quantum 
tunneling hence changing the local dipole moment of the 
vacuum, and the charges move in the medium due to the 
residual electric field giving rise to a conductive current. 
Following we name displacement current, jA, and 
matter current, jM- 


j = jM+jD- ( 10 ) 

Here Jm is a colored generalization of the usual electric 
current density which in a continuum notation is given 

by 


3 m= 9Y^Q3 c [ ^P^fjcip) ■ (11) 

..C J PO 

The displacement current arises from the polarization of 
the vacuum due to the decay of the electric field by the 
Schwinger mechanism: more precisely it is given by the 
time derivative of the local dipole moment induced by 
the particles pop-up, in the same way a time variation 
of the local dipole moment in a medium gives rise to 
a change in the local electric field [nl. According to 
the quantum tunneling interpretation of the Schwinger 
effect [3^ the dipole moment can be computed as 2 po /E 
where po corresponds to the kinetic energy of the particles 
coming out from the vacuum; taking into account Eq. (ED 
jo can be written, in the reference frame where particles 
are produced with vanishing longitudinal momentum, as 


JD = 


Nf 3 


EE 


d^p dNjc 2pT 
Po dr E ' 


( 12 ) 


where Nf corresponds to the number of flavors in the cal¬ 
culation. The color charge and current densities depend 
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on the particle distribution function: hence they link the 
Maxwell equations (I5|) to the kinetic equation ([51). We 
solve self-consistently the field and kinetic equations: in 
this way we take into account the back reaction of parti¬ 
cle production and propagation on the color field. 

At variance with the standard use of transport theory, 
in which one fixes a set of microscopic processes into the 
collision integral, we have developed an approach that 
fixes the total cross section in order to have the wanted 
77 /s. By means of this scheme we are able to use the 
Boltzmann equation to simulate the dynamical evolution 
of a fluid with specified shear viscosity, in analogy to 
what is done within hydrodynamical simulations. 

We use the CE approach to relate shear viscosity to 
temperature, cross section and density which is in agree¬ 
ment with Green-Kubo correlator results [^ . Therefore, 
we fix 77/5 and compute the pertinent total cross section 
by mean of the relation 

1 T 1 

<7tot = 3 j—r 7 — , (13) 

5 p 5 (a) 77 /s 


which is valid for a generic differential cross section 
da/dt ^ oi1/{t — as proved in [h^. In the above 

equation a = mo/'^.T, with tud the screening mass regu¬ 
lating the angular dependence of the cross section, while 


5(a) = 


{y'^+l)K3{‘2y)-yK2{2y) 


hii 


(14) 


where Kn is the Bessel function and the function h re¬ 
lates the transport cross section to the total one crtr{s) = 
atotHmh/s) being h(C) = 4C(l-kC) [(2C + 1 )^?t-( 1 + 1/C)- 
2]. The 5 (a) is the proper function accounting for the 
correct relaxation time = g{a)utotP associated to the 
shear viscosity transport coefficient. For isotropic cross 
section, i. e. mo —t 00 , the function g{a) is equal to 2/3 
and Eq. CSl) reduces to the relaxation time approxima¬ 
tion with T~^ = = <ytrP, while for finite value of ttid, 

which means anisotropic scatterings, 5 (a) < 2/3. 

We notice that, in the regime where viscous hydro¬ 
dynamics applies, the specific microscopic details of the 
cross section are irrelevant, and our approach is an effec¬ 
tive way to employ transport theory to simulate a fluid 
at a given g/s 1 ^. 


IV. FLUX TUBE DECAY IN A STATIC BOX 

In this Section we study the chromoelectric flux tube 
decay in a static box. We assume the box is a cube with 
side of 5 fm. Moreover we assume periodic boundary 
conditions for particles propagating in the box. This case 
is of academic interest, nevertheless it is useful because it 
allows us to introduce concepts which will be useful when 
we will consider more interesting case of a longitudinally 
expanding background, as well as it provides a further 
test of the numerical solution of Eq. ®. 


Given the symmetry of the problem, assuming at ini¬ 
tial time a homogeneous electric field along the z di¬ 
rection, then the system will remain homogeneous along 
the whole dynamical evolution: the electric field at later 
times as well as the currents and the invariant distribu¬ 
tion functions will depend only on time and not on space 
coordinates. Moreover it is easy to verify by Maxwell 
equations that neither magnetic fields nor transverse 
components of the electric field can develop during the 
evolution. Within these assumptions the evolution equa¬ 
tion for the classical field relevant for our problem is given 
by 


dE , ^ ^ 

n = -«'> ■ 

In our simulations in the static box case we use Eq. (CD 
to create particle pairs from the color electric field: at 
each time step, the value of E and of the volume box be¬ 
ing given, we compute the expected pair number, J\f, inte¬ 
grating Eq. CD over the volume box, then we distribute 
the M pairs uniformly in the box and with transverse 
momentum pT according to the distribution in Eq. CD; 
since the pairs have to pop out from the vacuum with 
vanishing total and longitudinal momenta, given px we 
extract randomly the azimuthal angle 4 > which uniquely 
determines Py = px sin (j) and Px = px cos (j) of one of the 
particles in the pair; finally the momentum direction of 
the second particle is given by tt — 0 . 

In Fig. [T]we plot the early times evolution of the color- 
electric field for the cases of calculations at fixed total 
cross section (upper panel) as well as fixed 77 /s (lower 
panel). In the upper panel of the figure the line labelled 
as ideal gas means a calculation with zero cross section. 

From the results in the upper panel of Fig. [1] we notice 
that for small coupling among the particles, the electric 
field rapidly decays through the Schwinger mechanism 
then it evolves with damped oscillations: the smaller the 
coupling is, the less efficient the damping is. Eventually 
for strong enough coupling the oscillations disappear and 
the color-electric fields just decays according to a power 
law. This dependence on the coupling strength is very 
easy to understand and is in last analysis due to the de¬ 
pendence of the electric conductivity of the plasma on the 
particle interaction strength [Wi| . In fact at initial times 
the particles are produced by the Schwinger effect with 
zero longitudinal momentum, then in case the coupling is 
small they are accelerated by the electric field thus gen¬ 
erating an electric current jm which at first gains energy 
from the field thus lowering its magnitude; because of 
the small coupling particles’s momenta is not random¬ 
ized thus currents are efficiently produced by the field. 
At some point of time the field is zero but the current 
Jm is still nonvanishing and positive hence causing a sign 
flip of the electric field and a negative acceleration of the 
charges, resulting eventually in Jm < 0 and an increase of 
the electric field. This process causes the field oscillations 
we observe in Fig. [TJ 

On the other hand, for large values of the coupling the 
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FIG. 1 : Early times evolution of the color-electric field for the 
cases of calculations at fixed total cross section (upper panel) 
as well as fixed ry/s (lower panel). 


scattering processes among the particles randomize mo¬ 
menta hence causing jm ~ 0 and the time evolution of 
the field is a pure decay by the Schwinger effect. To con¬ 
sider this point more closely we write down the analytical 
form of Equation (HU, namely 


JD 


C( 5 / 2 ) 


4-^2 

4 


3 Nf 3 


C—1 j — ^ C—1 


( 16 ) 


For example in the case one considers only the decay into 
gluons and assuming for simplicity Cg = 0 in Eq. one 
has for Jm = 0 the equation 

§ = , (17) 

where we have introduced the quantity 


ctd 


C(5/2)4-y2 

47r3 4 


E \9Q0e\^^^ 1 


C=1 


( 18 ) 


which is related to the Schwinger effect and is indepen¬ 
dent on the coupling among the particles. The solution 
of Eq. (IT 71 ) can be found as 


4 E( 0 ) 

(2 -I- 


( 19 ) 



FIG. 2: Time evolution of the ratio Pl/Pt for the static box, 
for several values of ri/s. 


showing that in absence of a conductive current the decay 
of the chromoelectric field is purely power law. 

In the lower panel of Fig. [T] we plot the early times evo¬ 
lution of the color electric field keeping 77/s fixed. The 
calculations have been performed assuming an isotropic 
cross section, and the relation among the total cross sec¬ 
tion and 77/s is computed by the CE method which, 
as proved in (^ . is in agreement with the relaxation 
time approximation for the case of an isotropic differ¬ 
ential cross section. From the qualitative point of view 
these results do not differ from those we obtain for calcu¬ 
lations at fixed cross section, and they can be understood 
in the same way the results we obtained in the case of a 
fixed cross section. 

In Fig. [ 2 ] we plot the time evolution of the ratio 
Pl/Pt for the case of the static box, for several val¬ 
ues of 77/s. As expected we find that for small values 
of 77/s the system is very efficient in removing the ini¬ 
tial anisotropy. We notice that the initial longitudinal 
pressure is negative and Pt/e = —1; indeed at initial 
time the system is made of pure color electric field for 
which T^’-' = diag (e, Pt, Pr, Pl) with e = P^/2, Pt = £ 
and Pl = —£■ On the other hand as soon as particles 
are produced by the Schwinger effect, they give a pos¬ 
itive contribution to the longitudinal pressure and the 
field magnitude decreases, eventually leading to a pos¬ 
itive pressure. Moreover because of the strong interac¬ 
tions among the particles produced, Pl and Pt tend to 
the same value within a very short time range, Tiso ~ 1 
fm/c. 

We notice that regardless of the value of 77/s the system 
tends to remove the initial pressure anisotropy, which is 
quite natural in the case of the static box since collisions 
always lead to the equilibrium state. However isotropiza- 
tion times are different and indeed in the case of a weakly 
coupled system, see for example the case Atttj/s = 10 in 
Fig. [21 the equilibrium state is reached in a much larger 
time than in the case of Airri/s = 1 . Moreover the os¬ 
cillations of the electric field strength in Fig. [1] lead to 
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many oscillations of Pl/Pt in the case ATTr]!s = 10, such 
oscilations being damped effectively within the first fm/c 
in the case A'Kr]ls = 1 . 


V. FLUX TUBE WITH BOOST-INVARIANT 
LONGITUDINAL EXPANSION 


In this Section we study the effect of a boost invariant 
longitudinal expansion on pair production from the de¬ 
cay of a chromoelectric flux tube. We assume the expan¬ 
sion takes place along the direction of the electric field; 
moreover the dynamics is invariant for boosts along the 
longitudinal direction. For this case we will discuss more 
results because it is closer to the description of the central 
rapidity region in the early stages of a relativistic heavy 
ion collision. Once again we closely follow the formalism 

of a. 

We remind very quickly the relevant equations in this 
case. The evolution equations for the color electric field 
are given by the pair of Maxwell equations 


dE _ dE 


( 20 ) 


where both current and charge densities are computed in 
the laboratory frame. Assuming boost invariance along 
the longitudinal direc tion imp lies that £z depends only 
on proper time r = \/t'^ — . We can combine the two 

Eqs. (Uni) to form a boost invariant equation, namely 


dE 

T— = zp-tj , 


( 21 ) 



FIG. 3: Chromoelectric field strength (main panel) and parti¬ 
cle number produced per unit of transverse area and rapidity 
(inset panel) as a function of time. The electric field is aver¬ 
aged in the central rapidity region |? 7 | < 0.5. 


To take into account the longitudinal expansion Eq. (HD 
has to be modifed as 

(23) 

where y denotes the momentum rapidity; TZ(jpt) depends 
only on transverse momentum and it is not affected for¬ 
mally by the expanding geometry; moreover we have in¬ 
troduced the two boost-invariant variables 


w = tpz - zpQ , V = Pot - zpz , (24) 


which can be rewritten as 
dE 

= psinhry — j coshry , ( 22 ) 

where the right hand side corresponds to the (minus) 
electric current computed in the reference frame where 
the time is the proper time, namely the local rest frame 
of the fluid. Equation (l2^ is in agreement with the boost 
invariant form of Maxwell equation used in , see also 

HI- 

To solve Eq. (1221) we adopt a finite difference scheme 
and prepare a box with a square cross section in the 
transverse direction, with —Xmax < x < Xmax and 
—Pmax < y < Vmax, and with cells in space-time rapidity, 
fixing the range of 77 in which we distribute the produced 
particles by -rjmax < V < Vmax with prnax = 2.5. This 
implementation corresponds to have a box with a longi¬ 
tudinal expansion since from the well known equations 
of relativistic kinematics one gets Zmax = ttanhrjrnax 
which, for Vmax sufficiently large, corresponds to a wall 
moving at ultrarelativistic speed along the longitudinal 
direction. We therefore distribute the pairs created by 
the decay of the flux tube with uniform probability in 
each of the cells in (x, y, y). 


and S(w) affects the longitudinal momentum distribution 
by forcing the condition 

z 

Pz = -po = Pt sinhy (25) 

for the produced pairs, with po = \/p\ + pi] such a con¬ 
dition is equivalent to assume that momentum rapidity 
of the produced pair is equal to space-time rapidity. The 
procedure we implement to create pairs in the case of 
the boost invariant longitudinal expansion is very simi¬ 
lar to the one we have described in the case of the static 
box, the only difference being that in the present case we 
iterate the static box procedure for each rapidity cell. 

A. Field decay, particle production and spectra 

In Fig. [3] we plot the color electric field strength aver¬ 
aged in the central rapidity region |y| < 0.5 (main panel) 
and particle number produced per unit of transverse area 
and rapidity (inset panel) as a function of time. The 
calculations are shown for three different values of y/s 
(calculations at fixed total cross section give similar re¬ 
sults); the relation among the total cross section and y/s 
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we have used in the simulation is the CE relation with 
an isotropic cross section, see Eqs. (ITITI) and (fHll . The 
electric field is averaged in the central rapidity region 
\y\ < 0.5 because in this region within a 10% of accuracy 
one has t ~ t with t corresponding to the laboratory 
frame in which the system expands, so a comparison with 
where the dynamics is followed in r is more meaning¬ 
ful. The initial value we have used in the simulations is 
E(t = 0) = 2.2 GeV^ in agreement with the LHC case of 
[54| but we obtain similar results by changing this value. 

As for the case of the static box we find that for the case 
of boost invariant longitudinal expansion the chromoelec¬ 
tric field experiences a rapid decay for small values oirj/s. 
Once again this is due to the fact that in this case the 
coupling among particles is large, meaning collisions are 
very effective in randomizing particle momenta in each 
cell hence damping conductive currents that might sus¬ 
tain the field. On the other hand for intermediate and 
large values of r;/s the electric field experiences stronger 
fluctuations during time evolution. 

In the inset of Fig. [3] we plot the number of produced 
gluons per unit of transverse area and rapidity versus 
time. We find that regardless of the value of rj/s we 
use in the simulation, the particles are produced at very 
early times, approximately within 0.5 fm/c, with the only 
expection of very few particles produced at later times in 
the case Airrijs = 10. We have checked that changing 
the initial value Eq of the electric field does not modify 
the production time in a considerable way unless Eq is 
very small, namely Eq 1 GeV^. Moreover the value of 
r]/s affects the conversion of the initial classical field to 
gluons only within a few percent: for example comparing 
the results for Airqls = 1 and Airyjs = 3 we find in the 
latter case a lowering of less than 10% on the number of 
particles produced. 

In the upper panel of Fig. 0] we plot the proper ki¬ 
netic energy density. Skin, at central rapidity (I 77 I < 0.5), 
versus laboratory time for three different values of 77 /s 
(calculations at fixed total cross section give similar re¬ 
sults). The relation among 77 /s and cross section is fixed 
by the CE relation with an isotropic cross section. We 
find that in the case of small 77 /s, which corresponds to 
the case of a strongly coupled system, the energy density 
decays asymptotically as Skin oc which is what is 

expected in the ideal hydrodynamic limit in the case of a 
one-dimensional expansion, in agreement with [d^. For 
the cases of larger 77 /s we find that a power law decay 
with a superimposed oscillation pattern is present. The 
thin dashed line in the figure corresponds to These 

results are in good agreement with those of [s^. 

In the lower panel of Fig. 0] we plot the plasma tem¬ 
perature as a function of time; it is obtained by data 
shown in the upper panel of Fig. 0] by assuming a perfect 
gas equation of state which gives T oc with pro¬ 
portionality constant being inversely proportional to the 
number of active degrees of freedom in the plasma. Our 
temperature is somehow larger than the one quoted in 
[ 5^1 because in the latter study both quarks and gluons 




FIG. 4: Proper kinetic energy density (upper panel) and local 
temperature (lower panel). All the quantities are averaged in 
the central rapidity region I?;! < 0.5 and viscosity has been 
fixed using an isotropic cross section. 


have been considered in the plasma, while in our case we 
only include gluons. The thin dashed line corresponds to 
which is the power law decay expected in the case 
of a one dimensional expansion of a non viscous fluid. 

In Fig. [5] we plot the gluon spectra at midrapidity 
\y\ < 0.5, for three different values of 77 /s. Upper panel 
corresponds to A-K-qls = 1, middle panel to Airrijs = 3 
and lower panel to A^rj/s = 10. For each value oirj/s the 
spectrum at three different times is shown. The thin solid 
black line corresponds to a thermal spectrum, namely 


dN 

prdprdy 


oc ptC ; 


(26) 


the above relation describes a thermalized system in three 
spatial dimensions at the temperature T = 1//3. In the 
figure the thermal spectrum is computed by taking the 
temperature at t = 1 fm/c from data plot in Fig. 0] . We 
find that for Airr]/s = 1 the system efficiently thermalizes: 
in fact the spectrum at t = 1 fm/c is of the form (1251) with 
temperature given (within a 2%) by the result in Fig. 01 
as it is evident by comparing the thermal spectrum (black 
line) with simulation data (dot-dashed thin green line). 
For Airp/s = 3 and AnTj/s = 10 the spectrum we obtain 






FIG. 5: Gluon spectra at midrapidity \y\ < 0.5, for three 
different values of 77 /s. Upper panel corresponds to 47rr;/s = 
1, middle panel to 4-Kr\/s = 3 and lower panel to An-qls = 10. 
For each value of 77 /s the spectrum at three different times 
is shown. Black thin solid line corresponds to the thermal 
spectrum in Eq. (I26II at t = 1 fm/c with temperature from 
data in Fig. |4l 


in the numerical calculation is in some disagreement with 
the thermal spectrum at the temperature plot in Fig. 0] 
at t = 1 fm/c, meaning the system is not completely 
thermalized in three dimensions. Moreover the very mild 
change in the slope of the spectrum we measure from 



FIG. 6 : Proper kinetic energy density (upper panel), lo¬ 
cal temperature (middle panel) and the ratio Pl/Pt (lower 
panel) as a function of time. All the quantities are averaged 
in the central rapidity region I 77 I < 0.5 and viscosity has been 
fixed using an isotropic cross section. 

t = 1 fm/c to t = 5 fm/c shows that the system does 
not cool down efficiently in this case, as it is expected 
because the large viscosity implies that a large part of 
energy dissipates into heat and the system cools down 
more slowly than the case of small viscosity. 

B. Pressure isotropization 

In this subsection we discuss pressure isotropization 
in the case with a longitudinal expanding geometry. In 
Fig. [5] we plot the ratio Pl/Pt where once again Pl and 
Pt correspond to the longitudinal and transverse pres¬ 
sure respectively. These quantities are computed cell by 
cell in the local rest frame of the fluid, then averaged 
in the rapidity range I 77 I < 0.5. The initial longitudi¬ 
nal pressure at initial time is negative and Pl/Pt = 
because at initial time the system is made of pure longi¬ 
tudinal chromoelectric field. On the other hand as soon 
as particles are produced, they give a positive contribu¬ 
tion to the longitudinal pressure and the field magnitude 
decreases, eventually leading to a positive pressure. For 
all the value of rj/s we consider in our simulations we find 
that the time needed to the total longitudinal pressure 
to be positive is about 0.2 fm/c. Moreover in the case 
47777/5 = 1 the strong interactions among the particles re¬ 
move the initial pressure anisotropy quite efficiently and 
quickly: in this case Pl/Pt = 0.7 within 0.6 fm/c, then 
the ratio tends to increase towards 1 within the time evo¬ 
lution of the system. This would justify the use of viscous 
hydrodynamics with tq ~ 0.6 fm/c as commonly done. 

On the other hand the larger the 77 /s of the fluid the 
larger the oscillations of Pl/Pt, compare for example 
the cases Air-rj/s = 1 and 47777/5 = 10 in Fig. [S] in the 
latter case Pl/Pt experiences several oscillations which 
follow the alternation of maxima of |i?| (corresponding 
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t [fm/c] 


FIG. 7: Ratios Pi/e, Pt/ s against time for s = 1 (upper 
panel), An-qls = 3 (middle panel) and Anrils = 10 (lower 
panel). All the quantities are averaged in the central rapidity 
region |? 7 | < 0.5 and viscosity has been fixed using an isotropic 
cross section. 


to minima of Pl since the field gives a negative contribu¬ 
tion to Pl) and zeros of E (corresponding to maxima of 
Pl)\ also, at large times Pl/Pt is quite smaller than 1. 
For intermediate values of 77 /s we still find some oscilla¬ 
tion which become less important for smaller viscosities. 
On the same footing the asymptotic value of Pl/Pt ap¬ 
proaches 1 when the viscosity becomes small. 


In Fig.[7]we plot the ratios Pl/s and Px/e asa, function 
of time. Here e corresponds to the total energy density, 
which takes into account energy density of both particles 
and field. We have shown results obtained for 47177/5 = 1 
(upper panel), Anrjf s = 3 (middle panel) and Anri/s = 10 
(lower panel). In Fig. [7] the red solid line corresponds to 
Pl/£ and the thin green line to Pt/£- The thin black 
dashed line corresponds to the conformal isotropic limit 
e = 3P. One can compare the results of Fig. [7] with 
those of 0 , where a classical Yang-Mills simulation with 
a 3 -|- ID expanding geometry is considered. The spirit of 
our comparison with Q is that, given we try to solve the 
same physical problem and we compute the same physi¬ 
cal quantities, a comparison of the final results obtained 
within the two approaches is interesting even if the theo¬ 
retical frameworks are different. In the weakest coupling 
case considered in Q, namely g = 0.1, Pl asymptoti¬ 
cally relaxes towards zero, which we might obtain if we 
introduce a larger viscosity than the one we consider in 
the present study. On the other hand the case g = 0.5 0 
produces Pl/s and Pt/s which lie in between our results 
for Airgls = 3 and Enrijs = 10. 


C. Anisotropic cross sections 

In this subsection we study the effect of changing the 
microscopic cross section from isotropic to anisotropic. In 
our collision integral we achieve this by tuning the Debye 
screening mass in the two body cross section, leaving the 
value of 77 /s fixed. In Fig. [5] we plot the time evolution 
of Pl/Pt for four different values of ttil) for the case of 
Attti/s = 1 (upper panel), Pkti/s = 3 (middle panel) and 
47777/5 = 10 (lower panel), and four different values of 
the Debye screening mass tttd. We remind that in our 
calculation the Debye mass is used as an infrared regu¬ 
lator of the differential cross section, and as a parameter 
which controls the anisotropy of the cross section: for 
very large values of ttid the differential cross section is 
isotropic, while for small values of mx) we get a forward 
peaked cross section. 

We find that lowering mu for a given value of 77/5 the 
plasma oscillations tend to be damped. This can be un¬ 
derstood because according to lowering mo while 
keeping fixed 77/5 amounts to increase isotropization of 
the distribution function; as a consequence the conduc¬ 
tive currents, which would sustain plasma oscillations in 
the late times evolution of the plasma, are damped. Nev¬ 
ertheless the effect on late time evolution of Pl/Pt is 
quite mild. 


D. Chapman-Enskog versus Relaxation Time 
Approximation 

In this subsection we study the effect of shifting from 
the CE to RTA when we relate 77/5 to the total cross 
section. It has to be noticed that our RTA does not cor- 
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FIG. 9: Ratio Pl / Pt versus time for ‘i'n’qjs = 3 and for mo = 
0.5 GeV (upper panel) and mo = 2.0 GeV (lower panel). In 
both panels we compare the Chapman-Enskog result (dashed 
red line) with the Relaxation time result (dot-dashed blue 
line). 


cross section [o^ . We limit ourselves to a particular value 
of ry/s, namely A'Kr]ls = 3 which is the intermediate case 
we have considered in the previous subsections. 


FIG. 8: Ratio Pl / Pt versus time for several values of m_D and 
Tj/s. All the quantities are averaged in the central rapidity 
region \r]\ < 0.5. Upper panel corresponds to Anri/s — 1, 
middle panel corresponds to Ai^rijs = 3 and lower panel to 
A^r]/s — 10 . 


respond to the RTA used in [H^: in fact in RTA 
corresponds to an ansatz for the collision integral in the 
Boltzmann equation; on the other hand in our calcula¬ 
tions, where we always solve the full collision integral in 
the Boltzmann equation, RTA refers only to an analytical 
equation which connects shear viscosity to microscopic 


In Fig. ini we plot Pl/Pt as a function of time for 
two values oi mu'- in the upper panel we consider rriD = 
0.5 GeV and in the lower panel tod = 2 GeV. In each 
panel, the red dashed line corresponds to the CE formula 
while the dot dashed blue line to the RTA formula. In 
the case of niD = 0.5 GeV the oscillations of Pl/Pt 
with CE are quite damped with respect to those of the 
RTA. This behaviour is easily understood by the results 
of [1^: as a matter of fact a given rj/s corresponds to a 
lower cross section for RTA, implying larger conductive 
currents which sustain plasma oscillations. In the case of 
Tod = 2.0 GeV the behaviour of Pl/Pt in RTA is similar 
to that in CE; in fact in this case the cross section is 
mostly isotropic, and it is known that RTA cross section 
is quite close to the CE one. 
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VI. SUMMARY AND CONCLUSIONS 

In this article we have reported about our results on 
the simulations of the Schwinger effect for color electric 
flux tubes, focusing on thermalization and isotropization 
of the fluid produced by the decay of the electric field in 
the tube. These simulations are important in the context 
of ultrarelativistic nuclear collisions, where flux tubes of 
strong color fields are expected to be produced in the 
early stage. According to the general understanding of 
high energy collision processes, this work is relevant both 
for heavy ion collisions, where a large number of tubes 
is expected, and for proton-nucleus and proton-proton 
collisions. 

We have studied both the cases of flux tube decay in 
a static box and in a box with a longitudinal expansion. 
The flux tube decay is described by the Schwinger ef¬ 
fect; the dynamical evolution of the field is assumed to 
be abelian so it satisfies Maxwell equations, and the dy¬ 
namics of the fluid produced by the decay is studied by 
the relativistic transport theory coupled to the Maxwell 
equations. We have formulated transport theory in terms 
of a fixed value of 77/5 rather than scattering amplitude; 
in this way we are able to describe the evolution of the 
fluid in terms of macroscopic quantities ( 77 /s) rather than 
insisting on a set of specific microscopic processes. The 
approach developed here, based on a stochastic solution 
of the relativistic Boltzman equation, has in perspective 
the advantage to be easily extended to 3 -I- ID realis¬ 
tic simulations, allowing to study the possible impact of 
the early dynamics on observables like radial and elliptic 
flow [^. Our study well matches where the same 
problem has been studied. The main difference between 
the present study and is that in the latter the relax¬ 
ation time approximation is used to simplify the collision 
integral in the Boltzmann equation, while in our study 
we do not use such an approximation hence solving the 
relativistic kinetic equation with the full collision inte¬ 
gral. 

We have considered the case of a system with boost in¬ 
variant longitudinal expansion, which is relevant for the 
case of simulation of early stages in heavy ion collisions. 
We have studied many aspects characterizing time evo¬ 
lution of the flux tube. In particular we have computed 
decay time of the color electric field, which we find to be 
a small fraction of fm/c; for large values of 77 /s we have 
found plasma oscillations develop within 1 fm/c and per¬ 
sist for several fm/c in agreement with [^; on the other 
hand for 77 /s < 0.3 plasma oscillations are very mild or 
even absent. 

We have also computed the time evolution of the lon¬ 
gitudinal to the transverse pressure ratio Pl/Pt- Given 
the particular initial field configuration and the absence 
of particle quanta at r = 0 we have Pl/Pt = —1 at 
initial time; however longitudinal pressure becomes pos¬ 
itive in r Ri 0.2 fm/c because particles are produced by 
the Schwinger mechanism, giving a positive contribution 
to the pressure. For 77 /s < 0.3 we have found isotropiza¬ 


tion times Tiso < 0.8 fm/c; for larger values of 77 /s plasma 
oscillations cause the system to be less isotropic. Our re¬ 
sults imply that for a fluid with A^rj/s < 3 isotropization 
time is less than 1 fm/c. This would justify the use of 
viscous hydrodynamics with initial times tq ~ 0.6 0.8 

fm/c for Attt]/ s = 1-^3. It remains to be understood if 
the oscillations of Pl/Pt for t < 1 fm/c affect some ob¬ 
servable in ultrarelativistic heavy ion collisions. Our ap¬ 
proach based on the stochastic solution of the Boltzmann 
equation is already set up to perform simulations start¬ 
ing from more realistic initial conditions and evolutions 
in 3-1-1 dimensions, which will be useful to understand 
quantitatively the effect of such early times behaviour of 
physical quantities on observables. 

We have also studied proper energy density and tem¬ 
perature evolution as well as particle spectra in the trans¬ 
verse plane for several values of 77 /s. We have found that 
in case of small 77 /s the dynamics is very efficient in both 
removing the initial anisotropy and in getting a quick 
thermalization; in fact for A-Krj/s = 1 we find thermal¬ 
ization time Ttherra < 1 fm/c by comparing the parti¬ 
cle spectrum at t = 1 fm/c with the thermal spectrum 
Eq. (1251) where we put the temperature equal to the one 
we obtain in the simulation, see lower panel of Fig. SI On 
the other hand, increasing 77 /s results in a mismatch be¬ 
tween the thermal and the actual spectra at t = I fm/c, 
meaning thermalization time is larger in these cases. 

To check the effect of microscopic details on the quan¬ 
tities we have computed, we have turned the microscopic 
cross section from isotropic to anisotropic by changing 
the value of the Debye screening mass in the two body 
cross section while keeping 77 /s fixed. We have found that 
changing mo does not lead to significant changes in the 
evolution of physical quantities, besides small deviations 
in the early times evolution and some mild damping of 
field oscillations in the case of the largest 77 /s. 

We have also considered the effect of turning from the 
CE to Relaxation Time Approximation (RTA) formula 
connecting the shear viscosity to the total cross section, 
still keeping the full collision integral in the Boltzmann 
equation. This is not exactly what has been done in 
where RTA is used as an ansatz for the collision integral. 
We have found that in the case of a forward peaked cross 
section changing from CE to RTA affects the behaviour 
of Pl/Pt, increasing oscillations at early times. Never¬ 
theless late times dynamics is almost unaffected by this 
change. 

There are several points we have not considered in this 
study. We have neglected rapidity fluctuations in the 
electric held which would break boost invariance. More¬ 
over in our description we have not included a longitu¬ 
dinal magnetic held at initial time which might modify 
sligthly the Schwinger amplitude [H^. Both huctuations 
in rapidity and the initial chromomagnetic held are quite 
interesting ingredients to add to the present study and we 
plan to add them in future works. Moreover it would be 
of a certain interest to consider geometries with several 
hux tubes pinned on the transverse plane, hence pushing 
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our calculations towards more realistic initial state condi¬ 
tions for heavy ion collisions as well as quantitative pre¬ 
dictions of observables of phenomenological interest like 
multiplicities, particle spectra and collective flows. We 
will consider these developments in forthcoming publica¬ 


tions. 
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